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Dynamic randomization and a 
randomization model for clinical 
trials data 

Lee D. Kaiser* f 

Randomization models are useful in supporting the validity of linear model analyses applied to data from a 
clinical trial that employed randomization via permuted blocks. Here, a randomization model for clinical trials 
data with arbitrary randomization methodology is developed, with treatment effect estimators and standard 
error estimators valid from a randomization perspective. A central limit theorem for the treatment effect 
estimator is also derived. As with permuted-blocks randomization, a typical linear model analysis provides 
results similar to the randomization model results when, roughly, unit effects display no pattern over time. A 
key requirement for the randomization inference is that the unconditional probability that any patient receives 
active treatment is constant across patients; when this probability condition is violated, the treatment effect 
estimator is biased from a randomization perspective. Most randomization methods for balanced, 1 to 1, treat- 
ment allocation satisfy this condition. However, many dynamic randomization methods for planned unbalanced 
treatment allocation, like 2 to 1, do not satisfy this constant probability condition, and these methods should be 
avoided. Copyright © 2012 John Wiley & Sons, Ltd. 
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1. Introduction 

A theory for estimation in a linear model, specifically accounting for the randomization used in a clinical 
trial will be developed. A particular goal is to compare the inferences from the usual population-based 
linear model with those from the randomization model. A further goal is to compare randomization 
model inference that follows from permuted-blocks randomization (PBR) with those from dynamic 
randomization methods. Here, 'dynamic randomization' refers to methodology in which the probability 
of assignment of a given patient to experimental treatment is a function of the patient's stratification 
variables and the stratification variables and treatment assignments of previously randomized patients. 
Other terminology has been used, such as adaptive randomization and minimization. 

Randomization of treatments to patients in clinical trials has perhaps most commonly been performed 
with PBR, often within strata defined by crossings of prognostic factors and/or investigational sites. 
Permuted blocks are used instead of 'complete randomization', which assigns the next patient to exper- 
imental treatment essentially by a fair coin toss, to ensure that the number of patients assigned to each 
treatment group remains similar within strata, regardless of the number of patients ultimately enrolled in 
the strata [1]. Efron [2] introduced biased-coin randomization, which assigns the next patient to exper- 
imental treatment with a biased coin if the counts of patients previously assigned to the two treatments 
are imbalanced. Pocock and Simon [3] generalized the method to include stratification variables, where 
the probability of assignment to experimental treatment depends on the balance of treatment counts 
across the margins of the stratification variables. There have subsequently been many proposals for 
other dynamic randomization procedures [4-7]. 
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A theory for estimation and testing exists for permutation-based randomizations paired with a linear 
model, where the linear model includes the randomization's stratification variables in addition to a treat- 
ment term, and leads to unbiased estimators of treatment effects and unbiased estimators of the variance 
of the treatment effect estimator [8-11]. The unbiasedness is over the distribution of the observations 
induced by repeated randomizations of treatments to fixed experimental units, under the assumption of 
unit- treatment additivity with fixed unit and treatment effects. 

There have been many theoretical developments for randomizations that are not of the permutation 
type, but these have mostly focused on randomization tests rather than estimation. 

o Smythe and Wei [12] derived the asymptotic distribution of a linear rank test of the null hypothesis 
of equality of two treatments when the assignment of treatments to patients is according to an 
urn-based randomization scheme. 

o Smith [13] summarized treatment allocation methods developed to date. For treatment allocation 
probabilities that are functions of (ri\ — «o)/(«o + ii)> where the «y's are the numbers of patients 
previously assigned to the two treatments, Smith developed an approximation to a randomization 
test based on the difference in sample means. Smith also developed similar approximations when 
treatment assignment probabilities are a particular function of stratification variables of the current 
and previous patients and developed asymptotic approximations to randomization tests. 

o Wei et al. [14] proposed a biased-coin scheme for the allocation of k treatments in proportion 
|i, . . ., tjfc, where %\ + ... + %k = \, and generalized the work of Smythe and Wei [12] to k >2. 

o Shao et al. [15] proposed a theory for clinical trials with dynamic randomization, but the results 
were developed on an assumption of a model that relates the dependent variable to the variables 
used in constructing hypothesis tests. 

o Rosenkranz [16] explored the sampling properties of a difference in sample means under 
various randomization schemes and noted that a Mest can be very conservative with biased- 
coin randomization. 

Regardless of the randomization method used, it is common in clinical trials data analysis to apply 
a statistical test that is based on assumptions of random sampling from a distribution. In general, these 
tests behave as follows for biased-coin and urn-based randomizations [2, 13, 16-18]. 

o If there is low-frequency variation in the responses of patients over time of entry into the trial, then 
a population-based test will tend to result in a larger /rvalue than the corresponding randomization 
test result. With this pattern of variation, patients enrolled closer together in time tend to respond 
more similarly than patients enrolled farther apart in time. 

o If there is high-frequency variation in responses over enrollment time, then a population-based test 
will tend to yield a smaller />-value than the corresponding randomization test result. In this case, 
patients enrolled closer together in time tend to respond less similarly than patients enrolled farther 
apart in time. This is generally noted to be less likely than the scenario above. 

o Finally, when responses appear to come from a homogeneous, time-independent model, then the 
two test results tend to be similar. 

These same qualitative comparisons apply for PBR when the blocking in time is ignored in the analysis 
[1], in which case the variation in patient responses is quantified by the intrablock correlation. 

Motivation for the present work on estimation and testing under a randomization model comes from 
a 'points to consider' document on adjustment for baseline covariates [19] that states '. . . techniques 
of dynamic allocation such as minimization are sometimes used to achieve balance across several 
factors simultaneously. Even if deterministic schemes are avoided, such methods remain highly con- 
troversial. Thus, applicants are strongly advised to avoid such methods.' Furthermore, Halpern and 
Brown [18] noted that dynamic randomization schemes 'can lead to complex randomization distribu- 
tions and hence to complex analyses about which little is currently known, and which, in fact, may 
present intractable difficulties.' With a randomization-based theory for estimation and testing applica- 
ble across all randomization methods, perhaps regulatory authorities will be more accepting of dynamic 
randomization methods. 

Section 2 contains the theory for randomization-based estimation in a linear model with unit-treatment 
additivity. Section 3 evaluates the estimation of the variance of the treatment effect estimator through 
a population-based linear model, in comparison with the randomization-based variance estimator. 
Section 4 compares the estimation results with randomization test results. Sections 5 through 7 explore 
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randomization-based inference when there is planned unbalanced treatment allocation, like two experi- 
mental to one control. A key requirement, often not met with proposed dynamic randomization methods 
for such unbalanced allocation, is that the probability that a patient receives experimental treatment 
should be constant across patients. The paper finishes with a discussion and conclusion. 

2. The randomization model 

Suppose n patients are in the trial, and they are assigned one of two treatments: experimental or control. 
The probability that a patient receives experimental treatment is a function of this patient's stratification 
variables and the stratification variables and treatment assignments of previous patients. No assumptions 
are made on this function just yet. Let 8 be the n x 1 vector that indicates assignment to experimental 
treatment: 8, = 0 if patient i receives control and 1 if patient i receives experimental treatment. Assume 
unit-treatment additivity, so that y, the n x 1 vector of observations on the patients, is 

y = yo + St, 

for unknown constant unit effects yo and a scalar treatment effect t. A vector of independent and 
identically distributed random variables could be added as an error term to the model, but this only 
complicates the development without changing the conclusions. 

Let X be an n x q matrix of covariates measured on the patients. Typically, some or all of the stratifi- 
cation variables used in the randomization will be included in X, but this is not required in what follows. 
Also, it will typically be thought that the variables in X are prognostic of patient outcome, but there is 
no assumed relationship between the unit effects in yo and the covariates in X. The only requirement is 
that 1, an n x 1 vector of ones, is in the column space of X and, without loss of generality, X has full 
column rank. As a simple example, X could have two columns, with the first column containing Is and 
the second column denoting Eastern Cooperative Oncology Group performance status of 1 or more. An 
analysis of covariance example has X with a column of 1 s and a column for a continuous covariate. 

We want to estimate x after adjustment for X, and so use the normal equations 

X'X X'8\ /p\ _ (X'y 
8'X &'8){x)-{8'y 

where P is a vector of nuisance parameters. With simple matrix algebra, 

8 'My 



x = 



S'MS 



, where M = I - X(X'X) _1 X' 



It is worthwhile to explicitly state what is fixed and what is random in this set-up, because this will 
provide the relative frequency basis for the inferences. The patients have entered the study in a given 
order, and while they could conceptually have entered in a different order, it is not clear what a reason- 
able probability model should be for these other potential patient entry orders. Also, would the patients' 
unit effects be different if they entered in a different order, and how should these effects be modeled? 
Would the patients' baseline covariates differ? Thus, patient entry order, the unknown unit effects, and 
patient covariates are considered fixed. What does clearly vary, and in a way that we can model, is the 
treatment assignment to patients, as captured in 8 . The variability in the observations y then comes only 
from the randomization mechanism, and the conceptual repetitions of the study arise from independent 
applications of the randomization method to these patients as entered with their fixed covariates and unit 
effects. This topic will be addressed further in the discussion section. 

Up to this point, no assumptions have been made on the randomization algorithm, and the only one 
needed to make progress is that E8 — pi, with p a constant: 1 tol randomization would have p — 0.5, 
2 to 1 randomization would have p — 2/3, and so on. The problems with randomization methods with 
E8 ^ pi are addressed in Section 6. 

This restriction to E8 — pi means that the unconditional probability that each patient receives 
experimental treatment is constant, even though the probability that a patient receives experimental 
treatment, conditional on the treatment allocations of previous patients, will often be different from p 
to more closely attain treatment balance, either overall or within stratification variable groupings. Note 
that it is not enough that E& averages to the desired fraction, like (l/n)l'£S = 1/2. With 1 to 1 treat- 
ment allocation, the E& — (1/2)1 condition is easily satisfied because of the usual symmetry present in 
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randomization algorithms, while dynamic randomization algorithms for unbalanced treatment alloca- 
tions often do not maintain a constant £8; value [20]. Examples of this are given in Section 6.1 and a 
potential solution for certain cases is given in Section 6. 

Note that t is a function of S'M8, which does not depend on unknown parameters and so is an ancillary 
statistic. In simple cases, as in Table I, S'MS is a function of sample sizes within the strata, and most 
would prefer to condition on observed sample sizes [21]. As shown in Appendix B, S'M8/« converges 
in probability to p(l — p) under mild conditions, so for large samples, inference that is conditional on 
S'MS and unconditional inference should be similar. However, for small samples, conditional inference 
may be preferred, so both approaches are developed in the next sections. 



2.1. Unconditional inference 

Appendix A.l develops approximations to the mean and variance of x under the restrictions of 
E& — pi and unit-treatment additivity and without conditioning on a given S'MS value. First, x is 
an approximately unbiased estimator of x. Second, 

„ ~ yo'MVMyp 

Var(x)% — (1) 

v ' (trMV) 2 

where tr() denotes the trace of a matrix and V = Var(8) is the variance-covariance matrix of 8. Var (x) 
is a function of the residuals in the regression of y 0 on X, so good use of covariance adjustment reduces 
the variance, which is as it should be. It is not so clear what the role of the randomization method is 
in the reduction of variance, because it is tied up in the V matrix. As an example, however, especially 
in a small study with a stratified analysis, a stratified permuted-blocks randomization is better than 
unstratified permuted-blocks randomization [22,23]. 
Next, it is shown in Appendix A. 1 that 

y'MVMy-x 2 S'MVM8 

Var(x) = ; ; (2) 

(trMV)S'MS - S'MVMS 

is an approximately unbiased estimator of the Var(x) approximation of Equation (1). Finally, under 
mild conditions it is shown in Appendix B that (t — i) / ^ Var (t) has asymptotically a standard normal 
distribution, so that tests and confidence intervals on x can be determined in the usual way. 

In a common class of randomization methods, t is unbiased, and not just approximately unbiased. 
Define a symmetrical 1 to 1 randomization method as one with P(S — S 0 ) = P(8 — 1 — S 0 ) for each 
treatment assignment sequence So - In words, the probability of obtaining any treatment sequence equals 
the probability of that sequence with the treatments switched between experimental and control. Most 
1 to 1 randomization methods treat the experimental and control arms symmetrically and satisfy this 
condition. For such randomization methods, the unbiasedness of t is shown in Appendix A.l. 

The diagonal elements of V are p{\ — p) and the off-diagonal elements are £(8,-8y) — p 2 , where 
£(8,8/) is the probability that patients i and j both receive experimental treatment. Analytical expres- 
sions for the first two moments of counts of assignments to treatments for urn-based randomizations 
[24] and for biased-coin randomizations [25] are available. For some cases, £(8,8 7 ) can be determined 
algorithmically by stepping through all assignment probabilities [26]. Alternatively, a straightforward 



Table I. S'MS for three examples of the X matrix. 



S'MS 



Intercept term only 

Two columns denoting membership in one of 
two strata, such as low ECOG status vs. high 
status 

Intercept and continuous covariate as in an 
ANCOVA 



(l/n 0 + l/"i) , 

where nj is the number of patients on treatment j 

(1/nio + l/'iii)" 1 + (l/n 2 o + 1/M21)" 1 , 

where rcy is the number of patients on treatment j in stratum i 

(l/no + l/nir^l-O/no + l/nir^xo-xOVS^y-x) 2 }, 
where nj is the number of patients on treatment j, Xjj is the 
covariate value of patient j on treatment i 
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solution for any randomization method is to rerandomize treatments to the patients in the trial according 
to the trial's algorithm millions of times and estimate V with the sample variance-covariance matrix of 
the resulting 8 s. 



2.2. Conditional inference 

For symmetrical 1 to 1 randomizations as defined above, inference conditional on a S'MS value is 
developed in Appendix B.l, where first it is shown that E (x|8'M8 = c) = x. Furthermore, letting 
Var c (t) = Var (t| S'MS = c), it is shown that 

yo'MV c Myo 

Var c (t) = , — 

' (S'MS) 2 

where V c = Var(S | S'MS = c). Finally, 

y'MV c My - x 2 trMV c MV c 



Var c (x) 



(S'MS) 2 - trMV c MV c 



is a conditionally unbiased estimator of Var c (t). 

Because most 1 to 1 randomizations are symmetrical, conditional inference is straightforward, as 
outlined above. For planned unbalanced allocation, however, even when ES — pi, Z?(8j | S'MS = c) is 
typically not constant over i. For example, consider PBR in blocks of 6 in a 2 to 1 ratio to 10 patients, 
withM = I-J/12. Unconditionally, £8; = 2/3. However, from Table I, if S'MS = (1/8+1/2)" 1 = 1.6, 
then £(8; | S'MS = 1.6) = 2/3 for 1 ^ i ^ 6, but E(8; | S'MS = 1.6) = 1 for 7 ^ i s= 10. 

The variation in the conditional probability of assignment to experimental treatment leads to 
conditional bias in x. Let p c = E(8 | S'MS = c), then 

8 ' My -'MS^Ux + ^ 



S'MS 

and the conditional bias is proportional to the covariance between p c and My 0 . The problem of con- 
ditional inference for planned unbalanced allocations with E8 — pi is illustrated in Section 5 and a 
potential solution is raised in the conclusion section. 



3. Comparison of standard error estimators — linear model versus 
randomization model 

Regardless of the randomization method used, it is common to analyze trial data with a population-based 
linear model 



y = Xp + Sx + e. (3) 

The elements of e are independent and identically distributed with zero mean and variance a 2 , 8 is con- 
sidered fixed, and X often includes stratification factors that are balanced by the randomization [27]. It 
is common, however, to ignore in the model the balancing in time that randomizations typically achieve. 

By construction, the treatment effect estimator is the same under the randomization model as under 
model (3); however, the variance estimator from Equation (2) will generally differ from the variance 
estimator derived from model (3). How similar are these estimators from a randomization perspective? 
The estimators will undoubtedly be different in any particular trial, but it is informative to consider their 
expectations. 



The linear model estimate of the variance of x is 



Var lm (x) = y' I-(X 8) (J* ' 



X'\ 1 



y'My 



(trM-l)S'MS trM-1 



8 j J (trM- 1)S'M8 

(4) 
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What is the expectation of Vari m (x) from a randomization perspective? Again, the numerator and 
denominator both contain random variables, but approximate by simply taking expectations individually 
of the numerator and denominator to obtain 

Z7,z~ ^ yo'Myo + x 2 trMV Vaiyc + x 2 

£VVari m (x) % 

' (trM-l)trMV trM-1 

yo'Myo Var r x 



(trM-l)trMV trM-1 

where the subscript r on E t and Var r denotes expectation and variance over the randomization method 
used. For what unit effects yo might Zs r Vari m (x) be similar to, less than, or greater than Var r (x)? Some 
guidance is given in the next two sections. 



3.1. Linear model versus randomization model standard error estimators when no pattern is expected 
for unit effects 

Suppose that apart from dependence on X, the unit effects yo display no further systematic variation. A 
reasonable model for this is yo = Xp + e, where e is a vector of independent and identically distributed 
random variables with expectation zero and finite variance a 2 . It is easy to see that the expectation of 
both E r Var lm (x) and Var r (x) equals a 2 /trMV. 

Therefore, if no discernable pattern in unit effects is expected, apart from variation owing to depen- 
dence on X, then the usual linear model estimate of the estimated treatment effect variance esti- 
mates a quantity similar to the proper randomization variance. This is true regardless of the type of 
randomization, as long as ES — pi. 

Note also that if complete randomization were used, with V = p{\ — p)\, then the randomization 
model variance estimator in (2) equals the linear model estimator in (4), illustrating the equivalence of 
randomization-based and population-based inference with complete randomization [28]. 



3.2. Linear model versus randomization model standard error estimators when unit effects have 
systematic variation 

The difference between the expectation of the linear model variance estimator in Equation (5) and the 
randomization variance in Equation (1) is proportional to yo'Myo /trM — yo'MVMyo/trMV. These 
quadratic forms have been evaluated in the context of accidental bias by Efron [2] and Smith [13] for 
biased-coin randomization. 

Efron defined accidental bias as the expectation, over the randomization used, of the squared bias of x 
calculated under an assumed regression model for the observations, where the regression model includes 
a covariate z not accounted for in the treatment effect estimator. This squared bias is proportional to 
z'MVMz, and therefore the work of these authors applies to the comparison of variance estimators con- 
sidered here. The applicable conclusion for biased-coin randomization is that when unit effect residuals, 
after regression on X, exhibit low-frequency variation or a smooth trend from patient to patient, then 
the linear model estimate of variance will tend to be larger than the randomization model estimate. With 
high-frequency variation, the comparison of the variances is reversed. 

Urn randomization is asymptotically free of accidental bias [17] so that for very large samples, the 
linear model variance estimate is likely to be similar to the randomization-based variance. For smaller 
samples, the control of sample size differences achieved by urn randomization is likely to lead to 
comparisons of variances similar to those for biased-coin randomization. 

This parallels the situation with PBR where the variance estimate ignoring the blocking in time will 
generally be greater than the estimate accounting for this blocking when there is low-frequency variation 
or a smooth time trend in the unit effect residuals, as quantified by the intrablock correlation [1]. 



4. Randomization tests 

The null hypothesis that treatment has absolutely no effect in any patient can be tested with a random- 
ization test. Under this hypothesis, the vector of observations y equals yo, and the observed x is just one 
of the possible treatment effect estimates that could be obtained by the random assignment of treatment 
labels to these patients. Is the observed x unusually large or small in this distribution of possible x values? 
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A randomization test is implemented by performing the randomization repeatedly with the patients 
enrolled in their fixed order and determining the resulting x r values. For an alternative hypothesis of pos- 
itive x, the randomization test />-value is estimated as the fraction of the x r values that are greater than 
the observed x; for a two-sided alternative, compare |x r | to |x|. As with estimation, the randomization 
test can be performed conditionally on 8'M8 = c or unconditionally. For the conditional randomization 
test, only those randomizations with 8 'MS = c, or suitably close to it, are used for the generation of the 
x r values. 

What are the first two moments of this randomization distribution of x r ? With E8 — pi, from 
Section 2.1, unconditionally E% % 0 and 

y'MVMy 

Var(x.) m — (6) 

V ; (trMV) 2 

Under the asymptotic normality conditions at the beginning of Appendix B, an approximation to the 
randomization test would be to compare t/yVar (x r ) with standard normal quantiles. This test statistic 
is similar to the Z-statistic obtained from use of SE (x) from Equation (2), except that a 'sum of squares 
for treatment' is removed in (2) relative to the variance in the randomization distribution. Through use 
of approximations as in Appendix A.l, the expectation of Equation (6) exceeds that of Equation (2) by 
approximately x 2 trMVMV/ (trMV) 2 , which is small for small x or large n. 

Thus, the randomization test inference will necessarily be similar to the inference that comes from 
simple use of x and its estimated standard error. As above, however, if the estimated standard error is 
obtained from a linear model based on X, then the randomization test need not agree with the linear 
model results. From Section 3.1, when the unit effects display no trend apart from dependence on X 
then the randomization inference should be similar to that from linear model inference. Likewise, from 
Section 3.2, when unit effects display a low-frequency trend, then the randomization test /7-value would 
be expected to be smaller than the />-value from the linear model, with the reverse order when the unit 
effects display high-frequency variation. This agrees with the summary of randomization tests given in 
the introduction. 



5. Conditional inference with planned unbalanced treatment allocation ratio 

Section 2.2 noted that planned unbalanced treatment allocation with E8 = pi yields an approximately 
unbiased treatment effect estimator unconditionally, but the estimator is typically biased conditioned on 
a 8'M8 value. This section illustrates by example the potential for bias with conditional inference and 
shows that the problem is applicable to PBR and to dynamic randomization. 

Consider the randomization of two treatments to 38 patients in a 2 to 1 ratio, with patients categorized 
by two binary stratification variables. The patient counts in the four cells defined by these two variables 
are respectively 10, 8, 11, and 9 for cells (0,0), (0,1), (1,0), and (1,1), with the first value of the pair 
indicating the Stratum 1 level and the second value indicating the Stratum 2 level. The order of entry of 
these patients to the simulated trial was determined at random and then held fixed. 

Two randomization methods are evaluated. The first is stratified PBR, with 2 to 1 allocation in blocks 
of six applied independently within each of the four cells. The second is a sequential method [5,29] that 
follows the recommendation at the end of Section 6. 1 to perform a symmetrical 1 to 1 to 1 randomization 
with three treatments, with two treatment arms subsequently combined (see Appendix C for details). The 
matrix X includes an intercept and main effects for the two strata. Each randomization method was repli- 
cated 10 million times to determine the V c matrices for this population. The methods were subsequently 
replicated one million times to determine the conditional properties of the estimators. 

For each randomization method, the expected value of the ancillary statistic 8'M8 is 8.3 and approx- 
imately 90% of the 8'M8 values are between 7.7 and 8.9; smaller values represent randomizations that 
generally exceed the 2 to 1 ratio and larger values represent randomizations with closer to 1 to 1 alloca- 
tion. Figure 1 illustrates the fluctuation across the patients in the empirically estimated E(&i |8'MS = c) 
for selected values of c, grouped to the nearest 0.2. For stratified PBR, £(8,- |8'M8 = c) is exactly 2/3 for 
earlier patients who belong to complete blocks, and similarly for the sequential method, E(&i |8'M8 = c) 
is approximately 2/3 for early patients. However, for patients in incomplete blocks with PBR and for the 
later patients with the sequential method, the conditional probability of assignment to active treatment 
diverges from 2/3, with £(8,- |8'M8 = c) generally above 2/3 for smaller values of c and below 2/3 for 
larger values of c. 
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Figure 1. E(8i|8 / M8 = c) versus patient randomization order, with values of c grouped to the nearest 0.2. 
(a) Stratified permuted-blocks randomization and (b) sequential randomization. See text for details on these 

randomization methods. 



From Section 2.2, the conditional bias is proportional to the covariance between p c = E(8\8'M8 — c) 
and the residuals of the unit effects after regression on X. Suppose here that unit effects are related to 
order of patient entry via a simple linear regression with a positive slope and coefficient of determi- 
nation R 2 , and let £ lm denote the expectation with respect to this model. Through simple algebra, the 
expectation of the conditional bias normalized by *J (£i m Var c (t)) is 



E lm E (x r | 8'M8 = c) - x 
J (£i m Var c (t)) 



p c 'Mz 



V{(trMV c )(l -R 2 )/R 2 + z'MV c Mz} 



where z is a unit-length vector proportional to patient randomization order. This normalized bias is dis- 
played in Figure 2 for each randomization method for R 2 — 0.5 and X = 1. The normalized bias is 
similar between randomization methods and is positive for smaller values of S'MS and negative for 
larger values of S'MS. At the smallest and largest values of S'MS displayed the normalized bias is 
approximately 10% in absolute value. 

Thus, for planned unbalanced treatment allocation with E8 — pi, randomization inference condi- 
tional on S'MS can be problematic, and this problem is present with both dynamic randomization and 
PBR. Because this problem arises in stratified PBR only for the patients in the last potentially incomplete 
block in each stratum, the bias is likely negligible in large trials. Likewise, for dynamic randomization, 
although based on limited empirical evaluations, p c appears to diverge from the desired fraction only for 
the latest patients randomized in a trial, so that the conditional bias for large trials would also seem to be 
negligible. 

The next sections shift to the problem of unconditional bias when Eh ^ pi and illustrate that this 
problem is common with planned unbalanced treatment allocation. 
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Figure 2. Conditional relative bias of x versus values of 8'M8 = c grouped to the nearest 0.2. Conditional relative 
bias is the conditional bias of x divided by the square root of the conditional variance of x. See text for details on 
the stratified permuted-blocks randomization and sequential randomization methods. 



6. Estimation and randomization tests when E8 ^ pi 

When ES — p ^ pi, then x is biased. Using a first-order approximation, the bias is 




S'Myp A _ £8'My 0 _ p'My 0 
8'MS ) % £8'M8 tr(MV) + p'Mp 



The bias is a function of the correlation between p and the residuals from the regression of y 0 on X. 

For the randomization test when ES ^ pi, the randomization distribution of the usual estimator x r is 
not centered at zero, but at p'My/{tr(MV) + p'Mp}, approximately. This nonzero expectation does not 
affect the Type I error rate by construction of the randomization test, although the power of the test may 
be affected [30]. 

A notable example of ES ^ pi is contained in the documents of the Endocrinologic and Metabolic 
Drugs Advisory Committee meeting held October 21, 2008 [30,31]. In this case, the linear model 
inference indicated a statistically significant treatment effect, while the randomization test did not, result- 
ing in p — 0.035 versus p — 0.06 for the analysis of one endpoint and p — 0.046 versus p = 0.15 for 
another endpoint. 

Because the bias is a function of p'My 0 , a simple method to obtain an approximately unbiased 
estimator of x is to augment X with the vector p, use [X p] in the normal equations, and proceed as before 
with estimation of x and the standard error. A randomization test would also proceed as before. The 
asymptotic normality results of Appendix B were developed with ES — pi, so it is uncertain whether 
asymptotic normality applies in this case. Even if it did, the best solution is not to use randomization 
methods with ES ^ pi in the first place. 



6.1. Dynamic randomization with planned unbalanced treatment allocation ratio 

Consider randomizations for a trial with 2 to 1 randomization of experimental to control treatment. 
/>(8,- = 1) = 2/3 for PBR applied in blocks of multiples of three. What about biased-coin randomiza- 
tion? A natural attempt is to use a biased coin whenever n e /2 differs too much from n c , where n e and n c 
denote the number of patients previously assigned to experimental and control. Suppose a threshold of 1 
or more is used, with P(8j = 1) = 0.9 for too many previous control arm assignments, i 5 (8 i = 1) = 0.2 
for too many previous experimental arm assignments, and P(8, = 1) = 2/3 otherwise. This procedure 
satisfies conditions CI and C2 in [13] and the probabilities satisfy the recent recommendation of Han 
et al. [6]. With the first patient's assignment to experimental treatment determined with probability 2/3, 
it is easy to calculate the P(8; = 1) values for i = 2, 3, 4, and 5 as 0.744, 0.467, 0.674, and 0.761. 
Stepping through the probabilities algorithmically yields P(82o = 1) = 0.775. The probabilities do not 
converge to 0.667, although the average of these 20 assignment probabilities is close to 0.667. 

The situation is worse with the marginal balance metric [6]. In this simple case, the biased coin is 
applied to the treatment that yields the smaller value of \(n e + l)/2 — n c \ versus \n e /2 — n c — 1|. 
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If the former is smaller, then experimental treatment is assigned with probability 0.9, and if the lat- 
ter is smaller, then experimental treatment is assigned with probability 0.2. With equal values and also 
for the first patient, experimental treatment is assigned with probability 2/3. P(8, = 1) values for i — 2, 
3, 4, and 5 are 0.433, 0.807, 0.881, and 0.342. Stepping through the probabilities algorithmically yields 
P(8 20 = 1) = 0.340. 

The generalization of 1 to 1 urn randomization to unbalanced allocation in Rosenberger and Lachin 
[32] is incorrect. For example, with 2 to 1 allocation, and starting with 2 E balls and 1 C ball, one 
application of their generalization would add 2 E balls if control treatment is assigned, with 1 C ball 
added if experimental treatment is assigned. The first patient is assigned experimental treatment with 
probability 2/3, and P(8,- = 1) values for i = 2, 3, 4, and 5 are 0.6, 0.590, 0.589, and 0.588. These prob- 
abilities converge to 0.586. A fix in this simple case is to add 4 E balls, rather than 2 E balls, if control 
treatment is assigned, with other aspects of the method unchanged. Then the P(8,- = 1) values for = 2, 
3, 4, and 5 are 0.620, 0.641, 0.656, and 0.662 and the probabilities converge to 2/3. In large samples, urn 
sampling behaves like complete randomization [17] so that deviation from the target allocation ratio is 
likely to be greater with urn sampling than with other methods. 

A solution for r\ : r 2 treatment allocation is to perform randomization for balanced allocation to 
r\ + r% treatment arms, and then combine t\ of these arms for the experimental arm assignment and ri 
of these arms for control. As long as the randomization method treats the n + ri arms symmetrically 
then E8 — pi will hold, and unconditional inference can be carried out as developed above. 



7. Discussion 

The estimation of x involves only 8, M, and y, and the estimator x is approximately unbiased as long as 
E& — pi and unit-treatment additivity holds. This approximate unbiasedness is true regardless of the 
relationship between y 0 and X and any other properties of the randomization algorithm. Therefore, at 
least from the point of view of unbiased estimation of x, there is no natural determination of X based on 
the randomization used, such as if the randomization is stratified on a particular factor, then that factor 
should be included in X. However, from the point of view of Var (x), covariates that are related to the 
unit effects should be included in X, as smaller residuals from a regression of yo on X reduce Var (x). It 
is unclear how manipulation of the randomization scheme influences Var (x), but it is likely that random- 
ization methods that seek to balance the covariates in X between treatment groups are somewhat better 
than those that do not, at least for small studies. Finally, for a proper randomization-based variance 
estimator, both X and the randomization method, through V, should be accounted for. On this point, the 
advice to 'analyze as you randomize' should be followed. 

Variability in the observations in the randomization model arises only from the random assignment 
of treatments to the patients in the order they entered the trial. This order and the patients' baseline 
covariates determine the distribution of the treatment assignment vector 8 and subsequent properties of 
the treatment effect estimator. In the discussion of randomization tests, this fixing of patient order is 
common [12, 13, 17, 26, 33]. A reviewer asked whether study conclusions can only be applied to these 
patients with the same order of entry? Now, it is common [34, 35] to consider that statistical inference 
in a randomized trial applies to the actual patients in the trial, with inference to the broader population 
of similar patients being either nonstatistical or based on assumptions that trial data follow specified 
distributions. In answer to the reviewer's question, it seems necessary to go a bit further and to say that 
the statistical inference from the randomization model applies only to these patients as they entered, 
at their specific stage of disease, and as treated at their specific study site. It is certainly plausible that 
treatment would have exerted a similar effect if patients had enrolled in a different order and under 
somewhat different circumstances, but to make statistical statements about such a treatment effect would 
require further assumptions. 

One generalization of the model of Section 2 is to drop the assumption of unit-treatment additivity. 
Let %i be the difference in response if patient i were given experimental treatment versus the patient's 
response with control treatment. The estimator t in Section 2 is in general a biased estimator of x, defined 
as the average of the x, values across the n patients. Using the first-order approximation to the expecta- 
tion of a ratio of random variables as in Section 2.1, the bias is a function of a contrast in the diagonal 
elements of the matrix VM. These diagonal elements will often be nearly constant, in which case the 
bias will be small. Gadbury [36] showed for permutation-based randomizations that the usual estimator 
of Var (x) is positively biased. If a similar result applied to the variance estimator in Equation (2) for the 
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more general randomizations considered here, then the randomization model estimators would yield, for 
example, nearly unbiased estimates of the average treatment effect, along with a conservative standard 
error estimate and conservative hypothesis tests on t. 

Another generalization is to more than two treatments. This appears to be relatively straightforward 
for a continuous endpoint with unit-treatment additivity through the introduction of a separate 8 vector 
for each of t — 1 treatments. 

It would be interesting to see whether randomization model inference could be extended from the 
continuous endpoints considered here to binary responses and time-to-event endpoints. Certainly, ran- 
domization tests of a null hypothesis of no treatment effect in any patient are available for such endpoints, 
but treatment effect estimation from a randomization perspective does not seem to be available. 

In addition to the analysis of the whole trial, it is important to explore treatment effects within sub- 
groups of patients. If E8 — pi and there is unit-treatment additivity within the subgroup, then a linear 
model estimate within the subgroup of this common treatment effect as in Section 2 is approximately 
unbiased, regardless of the covariates included in such an estimate. An estimate of the standard error of 
this estimator could come from Equation (2), with V determined only from this subset of patients. As in 
Section 3.1, if unit effects within this subset are expected to display no systematic variation other than 
dependence on the covariates, then the linear model estimate of the standard error should be a reasonable 
approximation to the randomization-based standard error estimate. 



8. Conclusion 

For symmetrical 1 to 1 randomization, inference conditional on S'MS and unconditional inference are 
both straightforward, and a randomization test will necessarily be similar to a test result obtained from x 
and its estimated variance, which is justified by a central limit theorem. With planned unbalanced treat- 
ment allocation, Proschan et al. [30] argued to 'minimize the use of minimization.' However, it is only 
with ES ^ pi that unconditional inference is problematic. Unfortunately, as illustrated in Section 6.1, 
many dynamic randomization methods for planned unbalanced allocation result in ES ^ pi, and care 
must be exercised to achieve E8 — pi, in which case unconditional inference is also straightforward. 
Although likely only relevant for smaller trials, with E8 — pi and planned unbalanced allocation, con- 
ditional inference will typically lead to bias in t. This conditional bias is also present with PBR when the 
blocking in time is ignored in the analysis. A potential repair to x is to include, as suggested in Section 6 
for unconditional inference, p c = Zs(8|8'M8 = c) as a covariate in X and proceed with inference based 
on V c = Var(S|S'M8 — c). This possible repair needs further evaluation. 

An additional similarity between PBR and dynamic randomization methods that attempt to bal- 
ance treatment counts is that inference from a population-based model will likely be similar to a 
randomization-based analysis when unit effects, or rather the residuals from the regression of unit effects 
on X, are essentially sampled from a homogenous population. Furthermore, for both PBR and dynamic 
methods, when unit effect residuals display low-frequency variation or a general trend related to enroll- 
ment time, then the estimated variance of x from a population-based model will tend to overestimate the 
randomization-based variance, with the opposite order in variance terms for high-frequency variation. 

In our company's experience, the use of dynamic randomization has never been an issue in disquali- 
fying a positive study. However, at the time of protocol review, we often receive comments from health 
authority review staff on our randomization method and trial analysis methods. It is very common that 
we receive a request for a randomization test for the main trial endpoints. I believe we would not receive 
such comments and requests if we routinely used PBR. With the results of this manuscript, owing to the 
similarity between randomization inference for dynamic randomization and PBR, it seems inconsistent 
to single out dynamic methods for the performance of a randomization test. Given that the blocking in 
time achieved by PBR and dynamic methods is routinely ignored in analyses, one should either (1) be 
satisfied that these analyses are likely to match a randomization analysis or to be conservative or (2) ask 
for a randomization test for any randomization method other than complete randomization. 

Appendix A.l 

Approximations to the unconditional moments of t will be developed here under the assumptions of 
E8 — pi and unit-treatment additivity. Recall that X is an n x q full column rank matrix with 1 in its 
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column space and that M = I - X(X'X) *X. First, note that 

„ = 8'M(y 0 + 8x) _ 8'Myp 
X_ 8'MS ~ 8'M8 

so that the moments of t can be found from the moments of 8'Myo/8'M8. Approximate this ratio with 
a first-degree Taylor series expansion about the means of the numerator and denominator to obtain the 
approximation 8'My 0 /trMV, where V = Var(8) and tr() denotes the trace of a matrix. Thus, 

8'My 0 \ _ £8'My 0 pl'My 0 



8'MS J trMV trMV 

because 1 is in the column space of X, so that I'M = 0'. Therefore, t is approximately unbiased. 
For the variance of t, 

Varl^Uvarf^ 
S'M8 j \ trMV J 

yo'MVMyo ... 

Var(x) , say. 



(trMV) 2 

For the estimation of Var (x), note that 

£y'MVMy = E (y 0 + St)' MVM (y 0 + 8t) 
= y'MVMyo + x 2 £8'MVMS 
= (trMV) 2 Var (t) + x 2 tr(MVMV) 



Define 



It is easy to see that 



y'MVMy - x 2 trMVMV 

Var 0 (x) = 

y ' (trMV) 2 - trMVMV 



£Var 0 (t) = Var (t) + |x 2 - (£x) 2 } 



trMVMV 



(trMV) 2 - trMVMV 



Because Ex % x, Var 0 (x) is an approximately unbiased estimator of Var (x). 

A problem with Varo (x), inherent in its unconditional nature, is that it averages over variation in 
8 'MS. A potentially better estimator that explicitly depends on 8 'MS can be found by requiring that 
the variance estimator match the linear model estimator of Equation (4) in Section 3.0 when complete 
randomization with V = p{\ — p)\ is used. With this requirement 

_ y'MVMVy - x 2 S'MVMS 
Var(T) ~ 8'MS trMV-S'MVMS 

Approximating the expectation of the ratio by the ratio of the expectations and approximating 
E (x 2 8'MVM8) by the product of the expectations, Var (x) is also an approximately unbiased estimator 
of Var (t). 

Appendix A. 2 shows for symmetrical 1 to 1 randomization methods, defined in Section 2.1, that t 
is unbiased conditional on 8'MS = c. The estimator t is therefore unconditionally unbiased for such 
randomization methods. 



Appendix A.2 

For the properties of t conditional on 8'MS = c, it is first necessary to evaluate E(& | 8'MS = c). 
In general, this conditional expectation is not proportional to the vector 1 even with E(8) — pi. An 
example for PBR in blocks of six for treatments in a 2 to 1 ratio is given in Section 2.2. 

In spite of potential bias issues with conditional inference with planned unbalanced treatment 
allocation, x is conditionally unbiased with symmetrical 1 to 1 randomization, defined in Section 2.1, 
where P(& = 8 0 ) = P(8 = 1 - 8 0 ). To show this, note first that if 8'MS = c then (1 - S)'M(1 - 8) = c. 
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Thus, let {8y} denote all treatment assignment sequences with Sy'MSy = c, where 8,2 = 1 — 8;i- 
Also, let /7,- ; denote the probability of obtaining the sequence 8,,. Because of symmetry, pn — pt2- 
By definition, 

E{8 | S'M8 = c) = ^ Pii 8ij = S ' ^'1^1 + (l-«.-i)> = (1/2)1 
v 1 ; P(8'M8 = c) PiS'ijMSij =c) K ' ' 

Hence, £ (t | 81M8 = c) = x for all c. 

Let Var c (x) denote Var (x | 8'M8 = c). Continuing with symmetrical 1 to 1 randomization, it is easy 
to see that 

y'MV c Myo 

where V c = Var(8 | 8'M8 = c). Following the approach in Appendix A.l, 

y'MVcMy - x 2 trMV c MV c 

Var c (x) = ; 

' (8'M8) 2 -trMV c MV c 

is a conditionally unbiased estimator of Var c (x). 



Appendix B 

Here, a proof is sketched of the asymptotic normality of x under the randomization model with E8 — pi. 
As in Section 2, x is adjusted for X, which has full column rank q and 1 in its column space. 

Intuitively, for reasonable randomization algorithms, treatment assignments for very large n have 
almost no dependence on treatment assignments for small n, so central limit theorems for sequences of 
random variables that satisfy a 'strong mixing' condition [37,38] should be helpful. However, the strong 
mixing conditions are difficult to verify for randomization algorithms. 

The approach taken here, then, is to approximate the sequence 8i, 82, ■ ■ ■ with an asymptotically 
equivalent triangular array {8„i, 8„2, • • •§««}• The original randomization scheme is applied within the 
array, except that treatment assignments are based on only the most recent m patients, with treatment 
assignments in earlier patients ignored. This makes the triangular array /^-dependent [38], and m will be 
allowed to increase with n . Note that permuted-blocks randomization generates a sequence of 8s that is 
already /w-dependent, where m is the block size. 

Because x = x + 8'Myo/8'M8 it suffices to show the asymptotic normality of ^/«8'Myo/8'M8. It 
will first be shown that 8'My 0 /^/n is asymptotically normal and then that 8'M8 jn ^ p p{\ — p), where 
— > p denotes convergence in probability. 

As before, let V = Var(8), with (i, j) element vtj . Make the following assumptions: 

(Al) yo'MVMyo//! -+ v > 0 as n ->- 00. 

(A2) Let (a„i , . . ., a nn )' = a„ = Myo and assume there is a finite A with |a„; | < A for all n and i . 
(A3) For all / , S y | Vjj | — > 5,- as n — > 00, and there is a B with Bj < B for all i . 

Let 8 M = (8„i, 8«2, • ■ -8«n)', the vector of the m-dependent array elements, where m is selected as the 
smallest integer less than or equal to w (°- 5_e ) for arbitrarily small e > 0. 

Now apply Berk's theorem [39]. His condition (i) on the uniform boundedness of E\a n i§ni\ 2+r 
applies for arbitrarily large r owing to Assumption A2. His condition (ii) on the uniform boundedness 
of Var(a„i8„j + . . . + a n j?> n j)/(j — i + 1) for 1 ^ i < j ^ n is satisfied because of both Assumptions 
A2 and A3. His condition (iii) on the limit of Var(a„i8„i + . . . + a nn ?> nn )/n follows from Assumption 
Al. Finally, his condition (iv) on lim„^oo m^ 2+2 ^ r ^ / n — 0 is satisfied for r > (1 — e)/e, where 6 was 
given above in the construction of the m-dependent array 8 n . It follows that (a M i8„i + . . - + a nn ?> nn )/y/n 
converges in distribution to A^(0, v). 

Next, assume 

(A4) the maximum eigenvalue of V, denoted A max , is o(n). 
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Using the inequality l'Vl ^ A max for any unit-length vector 1, it is then straightforward to show that 
8 'MS / n ^ p p ( 1 - p ) . First, let M — I - LL' , where L'L = \ q . Let 1, denote the i th column of L, where 
without loss of generality, li = 1/^/n. Then 

S'MS S'(I-J/h)8 £i>i(l^S) 2 



Now, 8'(I-J/n)8/« = 8'1/n - (8'1/n) 2 -> p p(l - p), because E(8'l/n) = 0 and Var(8'l/») = 
(l'Vl/»)/« ^ A max //7 — > 0. The result follows because for each i — 2, . . .,q, £(l'-8/ \/ri) — 0 and 
Var(l ■ 8 / y/n ) ^ A max /« — > 0. Thus, ^/« (x — x) converges in distribution to 7V(0, u/{/>(1 — p)} 2 ). 

Because v is unknown, practical application of the asymptotic normality of x for tests and confidence 
intervals requires a consistent estimator of v/{p(l — p)} 2 . It will now be shown that Var (x) of Equation 
(A1.3) provides one solution, but a further assumption is needed. 

(A5) A max , the maximum eigenvalue of V, is o(^/n). This is a stronger version of Assumption A4. 

Now, 

yo'MVMyo/rc - 2xy 0 'MVMS/n + (x 2 - x 2 ) 8'MVMS/n 



nVar (x) = 



(S'MS /«)(trMV/«) - S'MVMS /n 2 



The pieces of this expression will be handled in turn, starting with the denominator. As shown above, 
S'MS//! -+ p p(l - p). Next, with M = I - LL' as above, trMV/rt = (trV - E,-lJVL)/«. The trace 
of V is np(\ — p) and EjljVl;/« is bounded between 0 and qX max /n, which converges to zero by 
A5. Thus, (8'M8/n)(trMV/«) -+ p {p(l - p)} 2 . Continuing with the denominator, 8'MVMS/n 2 = 
(S'MVMS/S'MSjS'MS/n 2 ^ (A max />7)(8'M8/n) -> p 0 from both A5 and the convergence in 
probability of 8'M8/t7. 

For the numerator of /iVar(x), y 0 'MVMy 0 /« -* v by Assumption Al. Next, £(y 0 'MVM8) = 0 
and it will be shown that Var(y 0 'MVM8/«) 0, implying that y 0 'MVM8/« ->- p 0. The following 
repeatedly uses the inequality l'Vl ^ A max for any unit-length vector 1 and that x'Mx ^ x'x for any x. 

xz / ^tvM. / ^ yo'MVMVMVMyo y 0 'MVMVMy 0 ^ A max y 0 'MVMVMy 0 

Var(y 0 MVM8/n) = ^ 

VJ ' ; 72 2 y 0 'MVMVMy 0 n 2 



< A max yo'MV 2 M yo yo'MVMyp ^ A max 2 yo'MVM yo 
« 2 yo'MVM yo 



which converges to 0 because of Assumptions Al and A5. For the final term in the numerator of the 
expression, y/n (x 2 — x 2 ) converges in distribution to a central normal distribution, so in absolute value 



\Jn (x 2 -x 2 ) | S'MVMS S'MS 



|(x 2 -x 2 )S'MVM8//2| = 



S'MS n^/n 

^ y n x -x ) — 

which converges in probability to 0, completing the proof that (x — x) / ^Var (x) converges to a standard 
normal distribution. 



Appendix C 

The dynamic randomization method of Section 5 that randomizes patients in a 2 to 1 ratio proceeds as 
follows. Determine the stratification variable values si and 52 for the current patient to be randomized 
and determine «io, «n, and n 12 for previously randomized patients, where, for example, n \ \ equals the 
number of the previously randomized patients with Stratum 1 level equal to the level S\ of the current 
patient who were randomized to treatment level 1. Likewise, determine n 2 o, «2i> an d »22- The assign- 
ment probabilities for the current patient are determined by the balance in these counts first for Stratum 
1 and then for Stratum 2. If max(« ij) — min(n ij ) ^ 1 , then the treatment assignment probabilities depend 
on the pattern of then 10, «n, and n\2 counts as follows: 

- With three unique values, the probabilities are, respectively, 0.75, 0.15, and 0.10 for the treatments 
associated with the minimum, intermediate, and maximum treatment counts n \j 
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- With two values at the minimum count and one at the maximum, the probabilities are 0.40 for each 
of the treatments associated with the minimum count and 0.2 for the treatment associated with the 
maximum count 

- With one value at the minimum count and two values at the maximum, the probabilities are 0.12 
for each of the treatments associated with the maximum count and 0.76 for the treatment associated 
with the minimum count 

If max(«i/) — min(«i 7 ) = 0 then a similar check on Stratum 2 counts is performed, with unequal 
treatment assignment probabilities as above if there is imbalance. If there is balance of counts for Stratum 
2, then the treatment for the current patient is selected from 0, 1, and 2 with equal probability. 

To convert from 1 to 1 to 1 randomization to 2 to 1 randomization, patients randomized to treatment 

0 are given control treatment and patients randomized to treatments 1 or 2 are given experimental. 

The probabilities in the above algorithm were selected to match the conditional treatment assignment 
probabilities of the stratified PBR in blocks of six in the following sense. With complete randomization, 
the probability of assignment to experimental treatment is 2/3 regardless of previous treatment assign- 
ments, so the deviation from the target unconditional probability of 2/3 is zero. At the other extreme is 
a deterministic scheme where a fixed one-third of the patients are given control treatment and a fixed 
two-thirds are given experimental, so the squared deviation from the target unconditional probability of 
2/3 is (0 — 2/3) 2 /3 + 2(1 — 2/3) 2 /3 = 0.22. In the example of Section 5 the mean squared deviation 
as estimated from the simulation is approximately 0.05 for both the stratified PBR and the sequential 
method described above. 
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